<<<<<<< HEAD Introducción al análisis de secuencias NGS - Alineamiento

INTRODUCCIÓN

El paso de alineación o mapeo de secuencias cortas reads producto de la secuenciación NGS a un genoma de referencia se considera una parte integral de los análisis de genomas (Pham-Quoc, Kieu-Do, & Thinh, 2019).

El problema del mapeo de reads consiste en ubicar o alinear dentro de un genoma de referencia previamente secuenciado los millones de reads para luego ensamblarlas a un nuevo genoma (Bak, Bodziony, Migdałek, Pareek, & Żukowski§, 2020) o en el caso de que no exista un genoma de referencia realizar un ensamble de novo.

Esta carga de trabajo se la realiza con diferentes enfoques de alineación, sin embargo, al igual que cualquier software, cada enfoque presenta diferentes formas de ejecución y consecuentemente de análisis.

Para realizar esta tarea se han desarrollado varios programas como lo son, BWA (Li & Durbin, 2009), Bowtie 2 (JHU, 2016) y otros, cuya función principal es hacer coincidir secuencias cortas provenientes de NGS con secuencias de referencia (Bak, Bodziony, Migdałek, Pareek, & Żukowski§, 2020). Actualmente, BWA es un software muy popular para fines de alineamiento, debido a su compatibilidad con datos de Illumina y la integración de distintos algoritmos de alineamiento con mayor rendimiento para reads de distintas longitudes (Li & Durbin, 2009).

El resultado del mapeo o alineamiento se guarda en un archivo d etexto plano con formato SAM (Sequence Alignment Map) o su equivalente comprimido en binario (BAM). Este archivo contiene de forma codificada toda la información de dónde y cómo cada read se alineo al genoma referencia. De el se pueden luego extraer información de variantes para análisis poblacional o estudios de asociación genómica o información de conteo de reads en regiones codificantes para en análisis de expresión de genes.

FORMATO SAM

El formato SAM consta de un encabezado que comienza con el símbolo @ y una sección de alineamiento que contiene la información de cada uno de los reads que alineo al genoma de referencia. Para conocer más acerca del formato SAM puedes revisar el paper que describe el formato o simplmente en wikipedia. Los archivos SAM se pueden analizar y editar con el software SAMtools. A continuación se muestra un ejemplo del formato SAM con los 11 campos obligatorios del alineamiento y 1 campo opcional.

Formato SAM

Formato SAM

SOFTWARES DE BIOINFORMÁTICA

BWA

Es un algoritmo de alineación para alinear lecturas de secuencia o secuencias de consulta largas con un genoma de referencia grande como el humano. Elige automáticamente entre alineaciones locales y de extremo a extremo y admite lecturas de extremos emparejados (Li, 2013).

Puedes visualizar las opciones y más información de bwa Aquí

Samtools

Es un conjunto de utilidades que manipulan alineaciones en los formatos SAM (Sequence Alignment Map), BAM y CRAM. Convierte entre los formatos, clasifica, fusiona e indexa, y puede recuperar lecturas en cualquier región rápidamente, puedes tener más información sobre este programa Aquí

Objetivos del trabajo práctico

  • Realizar el alineamiento de una muestra en formato .fastq a un genoma de referencia.

Origen de las muestras

Para facilitar el proceso solo se utilizará una muestra la que provieme de la base de datos SRA del NCBI y corresponden a lecturas crudas del salmón del Atlántico Salmo salar en formato fastq, obtenidas por secuenciación de extremos emparejados con un secuenciador Illumina HiSeq2000.

La muestra se alineará al genoma de la mitocondria del salmón del Atlántico puesto que dado que su pequeño tamaño, el alineamiento no debería tardar más de unos pocos minutos.

Etapas del alineamiento.

  • Descargar Genoma de referencia de la mitocondria e indexar.

  • Alinear muestra con genoma de referencia.

  • Analizar el resultado del alineamiento.

Actividades

1 - Conectar a servidor Pomeo

Conecte al servidor pomeo usando la siguiente dirección IP: 200.54.220.141 puerto 22.

Para Windows puede conectarse con el software PuTTy y las claves de acceso dadas en la primera clase.

2 - Configurar bioconda e instalar programas para análisis

Para configurar el canal bioconda se debe ejecutar el siguiente comando

conda config --add channels bioconda

Instalar los software bwa

conda install -c bioconda bwa
    

Instalar samtools requiere un proceso de instalación un poco más largo, siga cada paso con atención.

conda install -c bioconda samtools

conda config --add channels bioconda 

conda config --add channels conda-forge 

conda install samtools==1.11

3 - Verificar directorios de instalación

Si desea constatar la ruta de instalación de un programa determinado puedes colocar el comando “whereis” en la terminal más el programa del que desea obtener la ruta.

 whereis sratoolkit 
 whereis samtools
 whereis bwa 

La salida de cada comando indicará la ruta de instalación como sigue, para cada programa del que desees obtener información.

/home2/usuario/miniconda3/bin/bwa
 

4 - Creación de directorio de trabajo y descarga de datos para alineamiento

Considerando que nuestras secuencias fastq originales tuvieron muy buena calidad, puedes trabajar directamente con ellas. Para ambos casos debes trasladar los archivos a una nueva carpeta, así:

Crear una carpeta denominada en tu usuario de home2

mkdir alineamiento

Ingresa a tu carpeta “alineamiento” y transfiere los archivos de clase anterior colocando los siguientes comandos en la terminal:

mv /home2/usuario/SRA_samples/SRR2006763/SRR2006763_1.fastq /home2/usuario/alineamiento/
mv /home2/usuario/SRA_samples/SRR2006763/SRR2006763_2.fastq /home2/usuario/alineamiento/

Lista tu carpeta de alineamiento para verificar que tienes lo necesario para el alineamiento, hasta ahora deben constar tus dos secuencias “SRR2006763_1.fastq” y “SRR2006763_2.fastq”

ls 

Ahora descargaremos el genoma de referencia de la mitocondria de Salmo salar en la misma carpeta de alineamiento, para ello:

En tu navegador entra al siguiente link

genoma_mitocondria

Encontrarás una tabla con el listado del genoma de referencia donde se incluyen todos los cromosomas y el genoma de la mitocondria y a cuyo identificardor RefSeq darás clic, como se ve en la siguiente imagen

Listado genoma Salmo salar

Listado genoma Salmo salar

Dar clic en la opción FASTA localizada bajo el título e identificador RefSeq de la referencia

Fasta option

Fasta option

Enviar la secuencia FASTA del genoma mitocondrial a un archivo

Genoma mitocondrial

Genoma mitocondrial

Ve al lugar de descargas de tu computador y encontrarás el archivo denominado “sequence.fasta” al que debes renombrar como “mt.fasta” y subirlo a POMEO, para ello puedes trabajar con WINSCP

Si trabajas con WINSCP debes iniciar sesión, guardar tus datos y conectarte como se indica en la siguiente imagen

WINSCP

WINSCP

Al conectar debes ingresar nuevamente tu contraseña y aparecerá la interfaz de tu servidor con las carpetas que tienes creadas, aquí ingresarás en tu carpeta de alineamiento y arrastrarás el archivo descagado del genoma mitocondrial a la misma.

Alternativamente puedes cargar el genoma iniciando sesión online en el servidor POMEO a traves del puerto 8787 precargado con Rstudio server.

http://200.54.220.141:8787/auth-sign-in

Cuando termines lo anterior puedes ingresar a POMEO, listar tu carpeta de alineamiento con ls y verás tu archivo “mt.fasta” junto con tus secuencias fastq.

5 - Indexación del genoma de referencia

Una vez que incluiste en tu carpeta de alineamiento todos los archivos descritos en pasos anteriores podemos proceder a la etapa inicial del alineamiento, que corresponde a la indexación del genoma de referencia con bwa usando el siguiente comando:

bwa index mt.fasta

La salida del comando dará como resultado 5 archivos con extensiones “amb”,“ann”,“bwt”,“pac” y “sa”

5 - Alineamiento

Para el alineamiento tendremos las siguientes etapas

  • Alineamiento de las secuencias contra el genoma de referencia, cuya salida será un archivo con extensión “.sam”

  • Conversión del archivo SAM a BAM

  • Inspeccionar el archivo .sam de salida

  • Ordenar lecturas alineadas por posición

  • Indexación con Samtools

  • Exploración de datos con Samtools

Para ejecutar todas las etapas anteriores en ese orden se debe crear un script con nano denominado aln_mt.sh

nano aln_mt.sh

En el script coloca las siguientes instrucciones

#!/bin/bash -l
# para alinear tus dos secuencias fastq al genoma mitocondrial
bwa mem mt.fasta SRR2006763_1.fastq SRR2006763_2.fastq > SRR2006763.sam 
# Transformar tu archivo sam a bam
samtools view -Sb -q 30 SRR2006763.sam > SRR2006763.bam 
# ordenar tu archivo binario bam
samtools sort SRR2006763.bam -o SRR2006763.sort.bam 
# indexar tu archivo bam
samtools index SRR2006763.sort.bam 

Ejecuta tu script con bash

bash aln_mt.sh

Ahora que ya tienes tus archivos SAM/BAM puedes observar tu archivo sam con el comando less de linux (recuerda que es un archivo de texto plano)

less SRR2006763.sam

también puedes realizar un análisis estadístico estandar con los siguientes comandos

samtools flagstat SRR2006763.bam > muestra_stat.txt

Investiga otros comandos de samtools para explorar el alineamiento.

Algunas observaciones adicionales

  • BWA-MEM se recomienda para secuencias de más de ~70 pb, ya que, generalmente BWA-MEM es más tolerante con los errores dadas las secuencias más largas.

  • Si quieres ver lo que sucede en cada etapa de alineamiento, puedes colocar en la terminal cada comando del script de forma individual, e.g.

    bwa mem mt.fasta SRR2006763_1.fastq SRR2006763_2.fastq > SRR2006763.sam # para alinear tus dos secuencias fastq al genoma mitocondrial samtools view -Sb -q 30 SRR2006763.sam > SRR2006763.bam # Transformar tu archivo sam a bam samtools sort SRR2006763.bam -o SRR2006763.sort.bam # ordenar tu archivo binario bam samtools index SRR2006763.sort.bam # indexar tu archivo bam

  • El archivo bam no se puede visualizar directamente porque está comprimido y en binario.

  • Procura listar tu carpeta cada vez que realices un trabajo, una descarga, o ejecutes un script para comprobar tus salidas.

Si quieres comprobar los tamaños de tus archivos al mismo tiempo que los listas, puedes hacerlo con

ls -l -h

REFERENCIAS

Pham-Quoc, C., Kieu-Do, B., & Thinh, T. N. (2019). A high-performance FPGA-based BWA-MEM DNA sequence alignment. Wiley, 1-12.

Li, H., & Durbin, R. (2009). Fast and accurate short read alignment with Burrows-Wheeler Transform. Bioinformatics, 1754-60.

JHU. (2016). Johns Hopkins University. Obtenido de Manual of Bowtie 2 - Fast and sensitive read alignment: http://bowtie-bio.sourceforge.net/bowtie2/manual.shtml

======= Introducción al análisis de secuencias NGS - Alineamiento

INTRODUCCIÓN

El paso de alineación o mapeo de secuencias cortas reads producto de la secuenciación NGS a un genoma de referencia se considera una parte integral de los análisis de genomas (Pham-Quoc, Kieu-Do, & Thinh, 2019).

El problema del mapeo de reads consiste en ubicar o alinear dentro de un genoma de referencia previamente secuenciado los millones de reads para luego ensamblarlas a un nuevo genoma (Bak, Bodziony, Migdałek, Pareek, & Żukowski§, 2020) o en el caso de que no exista un genoma de referencia realizar un ensamble de novo.

Esta carga de trabajo se la realiza con diferentes enfoques de alineación, sin embargo, al igual que cualquier software, cada enfoque presenta diferentes formas de ejecución y consecuentemente de análisis.

Para realizar esta tarea se han desarrollado varios programas como lo son, BWA (Li & Durbin, 2009), Bowtie 2 (JHU, 2016) y otros, cuya función principal es hacer coincidir secuencias cortas provenientes de NGS con secuencias de referencia (Bak, Bodziony, Migdałek, Pareek, & Żukowski§, 2020). Actualmente, BWA es un software muy popular para fines de alineamiento, debido a su compatibilidad con datos de Illumina y la integración de distintos algoritmos de alineamiento con mayor rendimiento para reads de distintas longitudes (Li & Durbin, 2009).

El resultado del mapeo o alineamiento se guarda en un archivo d etexto plano con formato SAM (Sequence Alignment Map) o su equivalente comprimido en binario (BAM). Este archivo contiene de forma codificada toda la información de dónde y cómo cada read se alineo al genoma referencia. De el se pueden luego extraer información de variantes para análisis poblacional o estudios de asociación genómica o información de conteo de reads en regiones codificantes para en análisis de expresión de genes.

FORMATO SAM

El formato SAM consta de un encabezado que comienza con el símbolo @ y una sección de alineamiento que contiene la información de cada uno de los reads que alineo al genoma de referencia. Para conocer más acerca del formato SAM puedes revisar el paper que describe el formato o simplmente en wikipedia. Los archivos SAM se pueden analizar y editar con el software SAMtools. A continuación se muestra un ejemplo del formato SAM con los 11 campos obligatorios del alineamiento y 1 campo opcional.

Formato SAM

SOFTWARES DE BIOINFORMÁTICA

BWA

Es un algoritmo de alineación para alinear lecturas de secuencia o secuencias de consulta largas con un genoma de referencia grande como el humano. Elige automáticamente entre alineaciones locales y de extremo a extremo y admite lecturas de extremos emparejados (Li, 2013).

Puedes visualizar las opciones y más información de bwa Aquí

Samtools

Es un conjunto de utilidades que manipulan alineaciones en los formatos SAM (Sequence Alignment Map), BAM y CRAM. Convierte entre los formatos, clasifica, fusiona e indexa, y puede recuperar lecturas en cualquier región rápidamente, puedes tener más información sobre este programa Aquí

Objetivos del trabajo práctico

  • Realizar el alineamiento de una muestra en formato .fastq a un genoma de referencia.

Origen de las muestras

Para facilitar el proceso solo se utilizará una muestra la que provieme de la base de datos SRA del NCBI y corresponden a lecturas crudas del salmón del Atlántico Salmo salar en formato fastq, obtenidas por secuenciación de extremos emparejados con un secuenciador Illumina HiSeq2000.

La muestra se alineará al genoma de la mitocondria del salmón del Atlántico puesto que dado que su pequeño tamaño, el alineamiento no debería tardar más de unos pocos minutos.

Etapas del alineamiento.

  • Descargar Genoma de referencia de la mitocondria e indexar.

  • Alinear muestra con genoma de referencia.

  • Analizar el resultado del alineamiento.

Actividades

1 - Conectar a servidor Pomeo

Conecte al servidor pomeo usando la siguiente dirección IP: 200.54.220.141 puerto 22.

Para Windows puede conectarse con el software PuTTy y las claves de acceso dadas en la primera clase.

2 - Configurar bioconda e instalar programas para análisis

Para configurar el canal bioconda se debe ejecutar el siguiente comando

conda config --add channels bioconda

Instalar los software bwa

conda install -c bioconda bwa
    

Instalar samtools requiere un proceso de instalación un poco más largo, siga cada paso con atención.

conda install -c bioconda samtools

conda config --add channels bioconda 

conda config --add channels conda-forge 

conda install samtools==1.11

3 - Verificar directorios de instalación

Si desea constatar la ruta de instalación de un programa determinado puedes colocar el comando “whereis” en la terminal más el programa del que desea obtener la ruta.

 whereis sratoolkit 
 whereis samtools
 whereis bwa 

La salida de cada comando indicará la ruta de instalación como sigue, para cada programa del que desees obtener información.

/home2/usuario/miniconda3/bin/bwa
 

4 - Creación de directorio de trabajo y descarga de datos para alineamiento

Considerando que nuestras secuencias fastq originales tuvieron muy buena calidad, puedes trabajar directamente con ellas. Para ambos casos debes trasladar los archivos a una nueva carpeta, así:

Crear una carpeta denominada en tu usuario de home2

mkdir alineamiento

Ingresa a tu carpeta “alineamiento” y transfiere los archivos de clase anterior colocando los siguientes comandos en la terminal:

mv /home2/usuario/SRA_samples/SRR2006763/SRR2006763_1.fastq /home2/usuario/alineamiento/
mv /home2/usuario/SRA_samples/SRR2006763/SRR2006763_2.fastq /home2/usuario/alineamiento/

Lista tu carpeta de alineamiento para verificar que tienes lo necesario para el alineamiento, hasta ahora deben constar tus dos secuencias “SRR2006763_1.fastq” y “SRR2006763_2.fastq”

ls 

Ahora descargaremos el genoma de referencia de la mitocondria de Salmo salar en la misma carpeta de alineamiento, para ello:

En tu navegador entra al siguiente link

genoma_mitocondria

Encontrarás una tabla con el listado del genoma de referencia donde se incluyen todos los cromosomas y el genoma de la mitocondria y a cuyo identificardor RefSeq darás clic, como se ve en la siguiente imagen

Listado genoma Salmo salar

Dar clic en la opción FASTA localizada bajo el título e identificador RefSeq de la referencia

Fasta option

Enviar la secuencia FASTA del genoma mitocondrial a un archivo

Genoma mitocondrial

Ve al lugar de descargas de tu computador y encontrarás el archivo denominado “sequence.fasta” al que debes renombrar como “mt.fasta” y subirlo a POMEO, para ello puedes trabajar con WINSCP

Si trabajas con WINSCP debes iniciar sesión, guardar tus datos y conectarte como se indica en la siguiente imagen

WINSCP

Al conectar debes ingresar nuevamente tu contraseña y aparecerá la interfaz de tu servidor con las carpetas que tienes creadas, aquí ingresarás en tu carpeta de alineamiento y arrastrarás el archivo descagado del genoma mitocondrial a la misma.

Alternativamente puedes cargar el genoma iniciando sesión online en el servidor POMEO a traves del puerto 8787 precargado con Rstudio server.

http://200.54.220.141:8787/auth-sign-in

Cuando termines lo anterior puedes ingresar a POMEO, listar tu carpeta de alineamiento con ls y verás tu archivo “mt.fasta” junto con tus secuencias fastq.

5 - Indexación del genoma de referencia

Una vez que incluiste en tu carpeta de alineamiento todos los archivos descritos en pasos anteriores podemos proceder a la etapa inicial del alineamiento, que corresponde a la indexación del genoma de referencia con bwa usando el siguiente comando:

bwa index mt.fasta

La salida del comando dará como resultado 5 archivos con extensiones “amb”,“ann”,“bwt”,“pac” y “sa”

6 - Alineamiento

Para el alineamiento tendremos las siguientes etapas

  • Alineamiento de las secuencias contra el genoma de referencia, cuya salida será un archivo con extensión “.sam”

  • Conversión del archivo SAM a BAM

  • Inspeccionar el archivo .sam de salida

  • Ordenar lecturas alineadas por posición

  • Indexación con Samtools

  • Exploración de datos con Samtools

Para ejecutar todas las etapas anteriores en ese orden se debe crear un script con nano denominado aln_mt.sh

nano aln_mt.sh

En el script coloca las siguientes instrucciones

#!/bin/bash -l
# para alinear tus dos secuencias fastq al genoma mitocondrial
bwa mem mt.fasta SRR2006763_1.fastq SRR2006763_2.fastq > SRR2006763.sam 
# Transformar tu archivo sam a bam
samtools view -Sb -q 30 SRR2006763.sam > SRR2006763.bam 
# ordenar tu archivo binario bam
samtools sort SRR2006763.bam -o SRR2006763.sort.bam 
# indexar tu archivo bam
samtools index SRR2006763.sort.bam 

Ejecuta tu script con bash

bash aln_mt.sh

7 - VISUALIZACIÓN DE ALINEAMIENTOS

Ahora que ya tienes tus archivos SAM/BAM puedes observar tu archivo sam con el comando less de linux (recuerda que es un archivo de texto plano)

less SRR2006763.sam

La salida de tu archivo sam se verá como este ejemplo:

Archivo .sam

  • Recuerda para salir del archivo puedes presionar la letra “q”

También puedes realizar análisis adicionales con los siguientes comandos:

Comando Descripción
samtools flagstat SRR2006763.bam > muestra_stat.txt análisis estadístico estándar
samtools flagstat SRR2006763.bam > muestra_stat.txt análisis estadístico estándar
samtools flags unmap proporciona los reads no mapeados
samtools flags 77 read 1 - emparejado no mapeado
samtools flags 141 read 2 - emparejado no mapeado
samtools flags 99 Reverse de un read 1 adecuadamente emparejado
samtools view -f 66 SRR2006763.bam | head -n 10 Busca solo reads emparejados en el archivo bam

Investiga otros comandos de samtools para explorar el alineamiento.

Algunas observaciones adicionales

  • BWA-MEM se recomienda para secuencias de más de ~70 pb, ya que, generalmente BWA-MEM es más tolerante con los errores dadas las secuencias más largas.

  • Si quieres ver lo que sucede en cada etapa de alineamiento, puedes colocar en la terminal cada comando del script de forma individual, e.g.

    bwa mem mt.fasta SRR2006763_1.fastq SRR2006763_2.fastq > SRR2006763.sam # para alinear tus dos secuencias fastq al genoma mitocondrial samtools view -Sb -q 30 SRR2006763.sam > SRR2006763.bam # Transformar tu archivo sam a bam samtools sort SRR2006763.bam -o SRR2006763.sort.bam # ordenar tu archivo binario bam samtools index SRR2006763.sort.bam # indexar tu archivo bam

  • El archivo bam no se puede visualizar directamente porque está comprimido y en binario.

  • Procura listar tu carpeta cada vez que realices un trabajo, una descarga, o ejecutes un script para comprobar tus salidas.

Si quieres comprobar los tamaños de tus archivos al mismo tiempo que los listas, puedes hacerlo con

ls -l -h

REFERENCIAS

Pham-Quoc, C., Kieu-Do, B., & Thinh, T. N. (2019). A high-performance FPGA-based BWA-MEM DNA sequence alignment. Wiley, 1-12.

Li, H., & Durbin, R. (2009). Fast and accurate short read alignment with Burrows-Wheeler Transform. Bioinformatics, 1754-60.

JHU. (2016). Johns Hopkins University. Obtenido de Manual of Bowtie 2 - Fast and sensitive read alignment: http://bowtie-bio.sourceforge.net/bowtie2/manual.shtml

>>>>>>> 9f5a27c820239e1c1dc27e82bc11b829dcd31bdc